clear all

file = "CleanedData/FHWACleanV4.csv";

xLabel = "yearRelativeEnd";
yLabel = "State";

data = readtable(file);
data = sortrows(data,xLabel);

data(data.SafetyClass~="Had Safety Inspections",:) = [];

%Clean out pre-1970 data
toDelete = isnan(data.LogRealMeanGasPrice);

data(toDelete,:) = [];

%Clean out repeated data points Registrations
toDelete = (data.StateCode == "CO" &(data.Year>=2002 & data.Year<=2009));
data(toDelete,:) = [];
toDelete = (data.StateCode == "CO" &(data.Year==2006));
data(toDelete,:) = [];
toDelete = (data.StateCode == "IN" &(data.Year==2006|data.Year==2007|data.Year==2009));
data(toDelete,:) = [];
toDelete = (data.StateCode == "MT" &(data.Year==2005));
data(toDelete,:) = [];
toDelete = (data.StateCode == "NJ" &(data.Year==2008));
data(toDelete,:) = [];
toDelete = (data.StateCode == "PR" &(data.Year==2002|data.Year==2003|data.Year==2005|data.Year==2008|data.Year==2009|data.Year==2010));
data(toDelete,:) = [];
toDelete = (data.StateCode == "TX" &(data.Year==2009));
data(toDelete,:) = [];
toDelete = (data.StateCode == "IL" &(data.Year==2011));
data(toDelete,:) = [];
toDelete = (data.StateCode == "NH" &(data.Year==2012));
data(toDelete,:) = [];
toDelete = (data.StateCode == "NY" &(data.Year==2012));
data(toDelete,:) = [];
%Clean out repeated data points gas use
toDelete = (data.StateCode == "RI" &(data.Year==2003|data.Year == 2005));
data(toDelete,:) = [];

min(data.Year)
max(data.Year)



xData = data{:,xLabel};
yData = data{:,yLabel};

figure('visible','off')
h = heatmap(data, xLabel, yLabel, 'ColorbarVisible', 'off');

XDisplayData = h.XDisplayData;
bins = [];

for i = 1:numel(XDisplayData)
    
    tempInd = str2double(XDisplayData{i});
    bins = [bins, tempInd];
    tempMod = mod(tempInd,5);
    if(~(tempMod==0))
        XDisplayData{i} = '';
    end
    
end
map = colormap([1,1,1;0,0,0]);
xLen = 20;
gr = (1+sqrt(5))/2;
yLen = xLen/gr;
f = figure('Units', 'centimeters', 'Position', [0,0,xLen,yLen], 'visible','off')
h = heatmap(data, xLabel, yLabel, 'ColorbarVisible', 'off', 'XDisplayLabels', XDisplayData, 'Colormap',map);
xlabel("Years Relative to Removal of Safety Inspections");
title("State vs. Year Relative to Event by Data Availability from 1970-2017 (Gas Use Model)")
savefig(f,"Plots/StateRelYearTwoWayGasUseModel.fig")
saveas(f,"Plots/StateRelYearTwoWayGasUseModel.png")


relYears = double(data.yearRelativeEnd);

counts = [];
yearCounts = [];
for year = bins
    tempCount = sum(relYears==year);
    counts = [counts tempCount];
    yearCounts = [yearCounts; year, tempCount, 16-tempCount]
end

baseVal = 0;

f = figure('Units', 'centimeters', 'Position', [0,0,xLen,yLen], 'visible','off')
area(bins, counts,baseVal, 'FaceColor', 'k')
hold on

linewidth = 2;

% plot([-6, -6], [0, 16], 'y', 'LineWidth', linewidth)
plot([-7, -7], [0, 16], 'r', 'LineWidth', linewidth)

% plot([-6, -6], [0, 16], 'y', 'LineWidth', linewidth)
plot([-7, -7], [0, 16], 'r', 'LineWidth', linewidth)

% plot([2, 2], [0, 16], 'y', 'LineWidth', linewidth)
plot([7, 7], [0, 16], 'r', 'LineWidth', linewidth)


legend("Number of States with Data", "Bound for 1 Missing State", 'Location', 'northwest')



xlabel("Years Relative to Removal of Safety Inspections");
ylabel("Number of States with Data");
title("Number of States With Data from 1970-2017 vs. Year Relative to Event (Gas Use Model)");
savefig(f,"Plots/NumStatesByRelYearGasUse.fig");
saveas(f,"Plots/NumStatesByRelYearGasUse.png");

dataAvailabilityTable = array2table(yearCounts, 'VariableNames', {'Years Relative To Treatment', 'Number of States with Data Available', 'Number of Missing States'});

writetable(dataAvailabilityTable, 'Tables/DataAvailabilityRelYearGasUseModel.csv')

